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Abstract 

DSTA  and  NTU  have  been  working  together  to  develop  a  methodology  to  simulate  the 
breakup  and  debris  throw  of  reinforced  concrete  (RC)  structures  under  internal 
explosion.  The  numerical  simulation  simulates  the  following  phases  of  responses: 
internal  explosive  loading;  structural  breakup  and  debris  launch;  debris  primary 
trajectory  and  post-impact  roll-and-bounce.  In  last  DDESB  seminar,  our  findings  on 
structural  breakup  process  were  reported  using  node-splitting  method  [1].  In  this  paper, 
the  focus  is  on  the  phases  following  the  structural  breakup.  As  a  result  of  a  structural 
breakup,  debris  pieces  of  different  sizes  are  formed  and  propelled  by  the  blast/gas 
pressures.  When  the  overpressure  subsides,  the  debris  pieces  fly  under  the  action  of 
gravity  and  the  air  drag.  Very  often,  the  debris  does  not  settle  at  its  first  landing 
position.  Instead,  it  impacts  the  ground  and  then  bounces  off  and  rolls,  or  even  breakup 
into  several  smaller  pieces.  This  paper  presents  the  investigation  on  all  these  post¬ 
breakup  processes,  including  the  effects  of  drag  coefficient  on  the  trajectory.  Results  are 
plotted  in  the  format  of  debris  dispersion  map.  Data  collected  from  the  Kasun  field  tests 
[2]  are  used  for  verification.  In  general,  the  simulation  results  are  in  good  accord  with 
the  field  test  results. 

1.  Introduction 

Debris  study  exists  in  many  fields.  Normally,  debris  study  includes  these  aspects:  firstly 
it  is  about  how  debris  is  produced  such  as  broken  part  from  space  shuttle  and  turbine, 
wind-borne  debris,  demolition  blast,  vehicles  transporting  explosives,  etc.  Secondly  it  is 
about  what  the  debris  flying  trajectory  is,  such  as  direct  prediction  of  the  recycling  site 
of  rocket  components  to  maximize  the  economic  benefits,  the  trajectory  of  debris  from 
demolition  blast  to  know  the  estimated  safe  clearance.  And  finally,  it  is  about  the  impact 
effect  of  debris  on  the  targets  normally  known  as  terminal  ballistics,  such  as  the  hit 
results  of  debris  from  shell  explosion,  the  hazards  of  turbine  hit  by  broken  debris,  space 
shuttle  survival  by  high  speed  space  debris,  wind-borne  debris  impact  on  building,  et  al. 

The  debris  involved  in  this  paper  is  from  the  break-up  of  reinforced  concrete  structure 
due  to  internal  explosion.  The  debris  trajectory  and  dispersion  are  studied  to  determine 
the  debris  hazard  zone  which  is  critical  in  determining  the  minimum  safe  clear  distances 
between  magazines  and  the  distances  between  inhabited  buildings  and  the  magazines 
concerned. 
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In  last  DDESB  seminar,  our  findings  on  structural  breakup  process  were  reported  using 
node-splitting  method.  Similar  methodology  for  the  simulation  of  masonry  wall  failure 
and  debris  scatter  due  to  high  explosive  loading  by  LS-DYNA  was  developed  [3].  In 
this  paper,  the  focus  is  on  the  phases  following  the  structural  breakup.  As  a  result  of  a 
structural  breakup,  debris  pieces  of  different  sizes  are  formed  and  propelled  by  the 
blast/gas  pressures.  When  the  overpressure  subsides,  the  debris  pieces  fly  under  the 
action  of  gravity  and  the  air  drag.  Very  often,  the  debris  does  not  settle  at  its  first 
landing  position.  Instead,  it  impacts  the  ground  and  then  bounces  off  and  rolls,  or  even 
breakup  into  several  smaller  pieces.  This  paper  presents  the  investigation  on  all  these 
post-breakup  processes,  including  the  effects  of  drag  coefficient  on  the  trajectory. 
Results  are  plotted  in  the  format  of  debris  dispersion  map.  Data  collected  from  the 
Kasun  field  tests  are  used  for  verification. 

2.  Debris  formation 

In  Quantity-Distance  (QD)  model  [4],  the  characteristics  of  the  debris,  like  the  initial 
velocity  and  mass  of  the  debris,  are  obtained  from  input  probability  distribution,  such 
as,  the  debris  mass  distribution  is  assigned  to  be  exponential.  The  debris  initial  velocity 
and  angle  distributions  are  assigned  to  be  normal.  While  in  present  study,  the  debris 
information  such  as  mass,  size,  ejection  position,  velocity  and  angle  are  retrieved 
according  to  certain  algorithms  (Seen  in  Appendix  A)  based  on  break-up  simulation 
output  by  LS-DYNA. 

3.  Debris  flying  trajectory 

3.1  Drag  force 

The  general  equation  for  the  drag  force  is[5] : 

F,=\c,PAV2  (3.1) 

where  Cd  is  drag  coefficient,  pa  is  air  density,  V  is  velocity  of  the  flying  object,  A  is  the 
object’s  cross-sectional  area.  In  this  study,  A  is  defined  as  debris’  maximum  sectional 
area  perpendicularity  to  the  motion  direction.  The  equation  means  that  the  drag  force  is 
proportional  to  the  drag  coefficient,  the  air  density,  the  normal  sectional  area,  and  the 
square  of  the  velocity. 

In  the  following,  a  new  parameter  k  is  introduced,  and  Eq.  3.1  is  re-written  as [6], 

Fd  =  mgkv2 ,  k  =  Cd  (3.2) 

2  mg 

3.2  Equations  of  debris  motion 

When  spinning  motion  is  excluded,  the  translational  motion  of  debris  can  be  described 
by  the  motion  of  its  centre.  The  forces  acting  on  the  debris  are  the  gravitational  force 
and  drag  force  only.  The  gravitational  force  acted  on  a  debris  mass  m  is  constant,  while 
the  drag  force  Fd  is  varying  and  is  proportional  to  the  square  of  the  velocity.  The 
kinematic  relationship  between  acceleration  a,  velocity  V  and  position  (x,y)  can  be 
expressed  in  differential  equations.  Hence,  the  debris  motion  is  governed  by  the 
differential  equations  as  follow  [6]: 
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dV  .  _  1tt2  dO  cos  0 
a  ~ ~~r  ~ ~g sm 9 - gkV  (p  =  —  =  -g—— 

dt  dt  V 


Vxm  —  =  Vcos0  V=^-  =  Vsin0 
dt  dt 


(3.3) 


Here,  V  is  the  velocity  of  the  debris,  0  is  the  slope  of  the  trajectory  to  the  horizontal,  g  is 
the  acceleration  due  to  gravity,  x  and  y  are  the  Cartesian  coordinates  of  the  debris.  The 
basic  parameters  involved  in  the  debris  trajectory  calculation  are  shown  in  Fig.  1: 


Fig.  1  Basic  parameters  of  trajectory  [6] 

Here,  Vo  and  do  are  the  initial  or  the  projective  velocity  and  angle  of  the  debris;  14  and 
Ok  are  the  landing  velocity  and  angle  when  debris  hits  the  ground;  Va  and  xa  are  the 
velocity  and  horizontal  displacement  when  debris  gets  to  the  vertex  of  the  trajectory;  H 
is  the  vertex  height;  L  is  the  debris  horizontal  displacement  before  it  hits  the  ground. 


To  solve  Eqs  (3.3),  we  divide  the  trajectory  path  into  a  series  of  small  intervals  [Oo,  0], 
then  the  whole  trajectory  could  be  traced  step  by  step  based  on  the  equations  3.4  [6], 


V(6)  = 


V0  cos  O0 


cos  djl  +  kvj  cos2  #0(/(#0)  -  /(#)) 


f(0)- 


sin# 
cos2  9 


.  i  On 
+  ln^|-  +  - 


(3.4) 


X  =  JC  +  Vp2  sin  2#0  -  V  sin  #  ^  +  V02  sin2  6>  -  V2  sin2  0_  ,  £  =  k(V02  sind0 +V2  S inO) 

0  2g(l  +  e)  ^  g(2  +  e) 

These  simple  analytical  formulae  of  Equations  3.4  enable  recurrent  computation  of  the 
main  parameters  of  motion:  H,  Va  ,  xa ,  L,  xa,  V*  ,  Ok ,  the  time  of  ascent  fa  and  the  total 
motion  time  t[6j. 

3.3  Determination  of  drag  coefficient  Q 

The  drag  coefficient  Q  is  a  number  that  aerodynamicists  use  to  model  all  of  the 
complex  dependencies  of  shape,  surface  roughness,  and  flow  conditions  on  moving 
object  drag.  This  dimensionless  value  quantifies  the  drag  or  resistance  of  an  object  in  a 
fluidic  environment  such  as  air  or  water.  The  higher  the  Q,  the  higher  aerodynamic  or 
hydrodynamic  drag  force  the  moving  body  will  experience  [7]. 

The  value  of  drag  coefficient  Q  is  a  function  of  the  shape  of  object  and  also  the 
Reynolds  number  (Re).  The  dimensionless  number  Re  is  defined  by  Re  =  plv/j)  where 
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l  represents  the  characteristic  length  scale  of  the  object  in  the  cross-sectional  plane  and  // 
the  dynamic  viscosity  of  the  medium.  In  the  present  study,  to  describe  the  debris  free 
flying  motion,  it  is  necessary  to  assess  its  range  of  Re.  The  air  density  p  is  assumed  to  be 
constant,  p= 1.293  kg/m3.  So  is  the  viscosity  of  air  (use  the  value  //=  1 .78x  1()  5  kg/(m-s) 
at  15°C  ).  As  such,  the  Reynolds  number  will  be  proportional  to  the  characteristic  size 
of  debris  (/)  and  the  velocity  (v).  To  estimate  the  lower  bound  of  Re,  set  /;  =  0.01m  and 
v i  =  10  m/s;  while  for  the  upper  bound,  set  I2  =  0.20  m  and  V2  =  500  m/s.  Thus  the  Re 
lies  in  the  range  of  [7.3xl03~7.3xl06]. 

Fig.  2a  shows  the  typical  relationship  between  the  drag  coefficient  Q  and  the  Reynolds 
number  Re  for  a  sphere  with  smooth  surface  [8].  As  can  seen  from  the  figure,  Q  is  less 
than  1  when  Re  is  greater  than  102;  and  C,j  is  about  0.5  when  Re  lies  in  the  interval  [103 
,  3xl05].  The  sudden  drop  of  Q,  when  Re  is  about  4xl05,  is  due  to  the  drag  crisis 
(which  is  associated  with  the  change  of  the  boundary  layer  on  the  sphere’s  surface  from 
laminar  to  turbulent).  From  [9]  it  is  observed  that  the  rough  surface  causes  the  drag 
crisis  shift  to  a  smaller  Re  interval,  accompanying  the  increase  of  the  minimum  Q. 


Fig.  2b  shows  curves  of  the  drag  coefficient  versus  Re  number  for  circular  cylinder  and 
square  cylinder  [9].  It  can  be  seem  that  the  square  cylinder  has  overall  drag  coefficient 
higher  than  the  circular  one.  Further,  the  square  cylinder  totally  skips  the  drag  crisis. 


In  the  present  debris  study,  the  formed  concrete  debris  will  be  in  multi-angular  irregular 
shape  with  serious  roughness.  To  adopt  a  reasonable  C d  value  for  present  study,  values 
for  three  benchmark  shapes  or  surfaces  are  further  considered  as  shown  in  Table  1. 

Table  1  Drag  coefficient  (’,/  [10] 
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The  values  for  the  spherical  shape  are  in  generally  the  smallest,  while  the  values  for  the 
circular  and  square  cylinders  are  higher.  As  the  debris  shape  will  never  be  a  perfect 
smooth  circular  cylinder  but  multi-angular,  it  makes  sense  to  assume  that  the  drag 
coefficient  for  debris  is  closer  to  that  for  the  square  cylinder,  whose  value  is 
independent  of  aspect  ratio  (L/D).  Again,  the  multi-angular  shape  of  debris  is  somehow 
between  the  circular  and  the  square.  Based  on  the  weakening  trend  of  drag  crisis  for 
rough  surface  and  totally  skip  of  the  drag  crisis  for  square  cylinder,  it  is  reasonable  to 
assume  that  the  drag  coefficient  for  the  debris  remains  more  or  less  constant  against  the 
Reynolds  number  in  the  our  interested  range  [7.3x10  -7.3x10  ].  Consider  the  value  Q 
=0.7  for  a  circular  cylinder  having  aspect  ratio  (L/D)  of  2,  and  the  value  Q  =2.0  for  a 
square  cylinder  of  any  aspect  ratio.  It  is  supposed  that  the  formed  debris  has  the 
transitional  shapes  from  the  smooth  circular  (Q  =0.7)  to  the  blunt  square  (Q  =2.0).  In 
the  absence  of  experimental  data,  values  of  drag  coefficient  for  the  intermediate  shapes 
are  assumed.  The  octagonal  shape  is  assumed  1.2;  while  the  hexagonal  is  1.6  as  shown 
in  [10].  So  the  range  (1.2  <  Q  <  2.0)  seems  fit  for  the  current  interest. 

3.4  Effects  of  Cci  on  debris  trajectory 

Fig.  3  shows  different  debris  trajectories  when  using  varied  drag  coefficient  Q  (0-2.0). 
The  debris  is  assumed  to  be  0.2  kg  having  launch  velocity  of  lOOm/s  at  an  angle  5°  with 
the  horizon.  It  can  be  seen  that  value  of  Q  has  significantly  effect  on  the  debris 
trajectory.  Lower  Q  leads  to  larger  travel  distance  and  it  is  more  prominent  when  Q  is 
less  than  1.0.  For  the  range  (1.2  <  Q  <  2.0)  of  current  interest,  the  travel  distance  when 
using  lower  one  (Q  =1.2)  is  about  101m;  while  using  the  higher  one  (Q  =2.0)  is  about 
81m,  which  is  approximately  81%  of  101m.  Empirically,  Q  =1.2  is  chosen  in  the 
current  study.  It  is  worth  noting  that  in  QD  model  the  input  for  the  drag  coefficient  is  a 
uniform  distribution  varying  between  1.0  and  2.0[4], 
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Fig.  3  Debris  Trajectories  by  various  (’,/  (m  =  0.2  kg,  V’,,  =  100  m/s,  d0  =  5°) 


3.4  Size  effect  on  C(i 

The  presence  of  drag  force  leads  to  size  effect  on  the  trajectory  because  mass  is 
proportional  to  volume  of  debris  while  drag  force  is  proportional  to  the  square  of  the 
square  of  the  cross-sectional  area.  With  the  same  initial  launching  conditions 
( v0  =  60 ml s ,  6{)  —  20° )  and  using  the  same  drag  coefficient  Cd  =1.2,  Fig.  4a  shows  the 

debris  travelling  distance  against  debris  mass.  The  relationship  appears  to  be  hyperbolic. 
For  smaller  debris,  the  travelling  distance  increases  almost  linearly  with  the  debris  mass, 
but  for  larger  debris,  the  rate  of  increases  is  seen  reduced  drastically  due  to  stronger 
drag.  Note  that  if  drag  is  absent  (Q  =0),  the  travelling  distance  should  be  the  same  for 
all  debris. 
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Fig.  4  Size  effect  on  Cd 
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(a)  Travelling  distance  vs.  debris’  mass 


To  illustrate  the  size  effect,  two  pieces  of  prototype  debris  are  used.  The  first  is  a  17-kg 
debris  and  the  second  one  is  a  (17/8)-kg  debris.  The  travelling  distances  for  each 
prototype  debris  and  its  half-scale  debris  for  a  series  of  launch  velocities  are  computed. 
Their  relationship  are  plotted  in  Fig.  4b.  Two  distinct  curves  are  seen.  It  illustrates  the 
size  effect  that  the  relationship  depends  on  the  mass  of  the  protype  debris. 

3.6  Debris  bounce 
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If  debris  hit  the  ground  at  a  shallow  angle,  it  will  bounce.  In  this  study,  debris  bounce  is 
considered  according  to  the  criterion  [11]  below. 


Fig.  5  Parameters  in  debris  bounce 

Vf  =  0 ;  ef  =  0  for  0.  >  13°  and  V(  >18 m/s  (3.5) 

Vf  =(1.0-  0.004760,  )Vt ;  Gf  =  0.4840,  otherwise 

Vt  is  the  velocity  of  the  debris  before  impact,  Qi  is  the  impact  angle,  measured  from  the 
ground  up,  Vf  is  the  post-impact  rebound  velocity  and  6f  is  the  post-impact  rebound 
angle  also  measured  from  the  ground  up. 

In  this  study,  the  total  travelling  distance  of  debris  is  limited  to  a  maximum  of  two 
bounces  (if  applicable)  since  the  effect  of  the  further  bounces  becomes  insignificant. 

3.7  Dispersion  map 

Debris  dispersion  map  can  provide  the  detailed  distribution  of  debris  settled  down  on 
ground.  It  reveals  the  debris  distribution  in  radial  and  circumferential  directions.  The 
dispersion  maps  are  presented  in  two  forms.  The  first  is  plotted  against  the  absolute 
debris  number.  The  second  is  plotted  against  the  normalized  number  (which  is  defined 
as  debris  number  per  m2  in  a  cell,  namely  debris  density).  The  debris  density  is  often 
used  as  a  parameter  to  evaluate  the  debris  hazard.  The  dispersion  map  can  be  used  to 
assess  the  hazard  zone  and  Inhabited  Building  Distance  (IBD). 

3.8  Hazard  zone  and  IBD 

The  ultimate  objective  of  this  study  on  debris  is  to  address  the  operational  need  with 
regard  to  the  safety/hazard  zone  around  a  magazine.  It  is  to  establish  some  guidelines 
for  identifying  safe  zone  against  the  hazardous  debris.  A  piece  of  debris  having  kinetic 
energy  (=  l/i  mv2)  of  more  than  79  J  is  considered  hazardous  or  lethal.  The  hazard  zone 
is  defined  as  a  patch  of  area  on  plan  such  that  the  lethal  threat  level  exceeds  one 
hazardous  piece  of  debris  per  56  m2[12],  A  hazard  zone  can  be  of  any  shape  as  long  as  it 
envelopes  all  hazardous  patches.  The  IBD  is  an  index  radius  of  the  smallest  circle  that 
encompasses  the  hazard  zone  totally  as  shown  in  Fig.  6. 
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Fig.  6  Hazard  zone  and  IBD 


4.  Case  study 


This  section  presents  the  illustrative  modeling  of  debris  throw  after  an  RC  magazine 
break-up  subjected  to  internal  blast  load.  Fig.  7a  is  the  filed  test  physical  object  in 
Kasun  II  tests  No.  2[2],  Fig.  7b  is  the  simulation  model  for  half  of  the  test  object. 
Internal  dimensions  of  the  cubic  box-shape  RC  structure  are  2.0m  x  2.0m  x  2.0m.  Wall 
and  roof  are  of  0.15m  thick  and  are  reinforced  with  0  12  mm  bars  at  100  mm  centers  in 
two  directions  having  double  layers.  The  charge  weight  is  20  kg-equivalent  TNT  having 
a  loading  density  of  2.5  kg/m3. 


4.1  Simulation  model 


Fig.7  Field  test  object  and  numerical  model 

This  physical  domain  is  modeled  by  finite  elements.  The  finite  element  model  includes 
discretization  of  the  solid  magazine  structure  and  the  air  space  inside  and  outside.  Due 
to  symmetry,  only  half- structure  model  is  built  to  save  calculation  cost  and  the 
symmetry  plane  locates  at  the  X-Z  plane  as  can  be  seen  in  Fig.  7b. 

4.2  Comparison  of  simulation  results  with  field  record  from  Kasun  Test 

4.2. 1  Comparison  of  launch  velocity  and  angle  distribution 

Launch  velocity  is  an  important  parameter  in  debris  hazard  study.  Figs.  8a  and  b  show 
the  launch  velocity  and  angle  distribution  respectively  obtained  from  simulation. 
Corresponding  results  obtained  from  field  test  are  plotted  in  Figs.  8c  and  8d.  Launch 
angle  is  measured  from  the  horizontal  plane.  Upward  is  positive;  downward  is  negative. 
The  launch  angle  varies  within  the  range  [-5°,  90°]  and  sorted  at  5°  intervals.  Note  that 
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if  drag  force  is  absent,  the  maximum  horizontal  travelling  distance  is  derived  from 
debris  having  launch  angle  of  45  degree.  With  the  presence  of  drag  force,  the  maximum 
distance  should  be  derived  from  debris  having  launch  angle  lower  than  45  degree. 


(a)  Launch  velocity  by  simulation  (b)  Launch  angle  by  simulation 


(c)  Launch  velocity  by  field  test  [13]  (d)  Launch  angle  by  field  test  [13] 

Fig.  8  Comparison  of  launch  velocities  and  angles 

It  can  be  seen  from  Fig.  8a  that  majority  of  debris  has  launch  velocity  in  the  range  of 
[50,  100]  m/s,  and  has  the  highest  density  (about  45%)  lies  in  the  range  of  [60,  70]  m/s; 
also  from  Fig.  8b  that  launch  angles  have  two  clouds  (71%  and  17%)  in  the  ranges  [0°, 
20°]  and  [80°,  90°],  respectively,  with  the  highest  density  (51%)  occurred  in  the  range 
of  [5°,  10°].  Further  scrutinizing  reveals  that  debris  having  low  launch  angle  are  mostly 
derived  from  the  walls  and  debris  having  high  launch  angle  are  mostly  derived  from  the 
roof.  Note  that  Figs.  8c  and  8d  are  the  probability  density  function  of  launch  velocity 
and  angle  deduced  from  high  speed  video  record  of  the  field  test,  but  limited  to  visible 
debris  having  launch  angle  not  more  than  45°.  It  can  be  seen  from  Fig.  8c  that  the 
velocities  are  nearly  of  normal  distribution  with  a  mean  velocity  of  about  80m/s,  which 
is  comparable  with  the  actual  [13]  mean  velocity  of  61m/s.  Comparisons  of  launch 
velocity  distributions  in  Fig.  8  show  that  simulation  results  agree  well  with  field  test 
record.  However,  the  launch  angle  distributions  only  poorly  agree  though  the  trends 
seemingly  match.  Nevertheless,  when  field  data  are  poorly  fit  to  a  normal  distribution,  it 
yields  a  mean  angle  of  7°,  which  agrees  well  with  the  simulation  results  and  also  the 
Klotz  group  recommendation  [13]  of  5°  with  a  standard  deviation  of  ±6°  (as  adopted  in 
the  prediction  tool,  namely  ‘KG-ET’),  having  95%  of  debris  in  the  range  of  [-7°,  +17°]. 

4.2.2  Comparison  of  dispersion  and  hazard  zone 


-9- 


34th  DoD  Explosive  Safety  Seminar,  Portland,  Oregon, 

13-15  July  2010 

Numbers  of  debris  pieces  finally  settled  on  ground  are  counted  in  19  equal  sectors, 
namely  Sector  A  to  Sector  S  over  a  190-degree  plan  view  as  shown  in  Fig.  9a.  The 
angle  of  each  sector  is  10  degrees,  with  the  first  one  starting  from  -5°  and  the  last  one 
ending  at  185°.  Along  the  radial  direction,  each  sector  has  10m  interval.  Fig.  9a  depicts 
the  collection  sectors.  In  the  real  field  test,  it  is  a  very  tedious  and  hard  work  to  collect 
all  the  debris.  To  alleviate  the  hardship,  the  debris  collation  is  not  exhaustive  and  only 
those  debris  pieces  having  mass  greater  than  0.055kg  were  collected  during  the  Kasun 
field  test.  Fig.  9b  plots  the  dispersion  map  from  the  Kasun  test  data.  Fig.  9c  shows  the 
debris  dispersion  map  obtained  from  simulation  which  includes  the  effect  of  bouncing 
twice  after  first  impact.  It  should  be  noted  that  the  smallest  minute-sized  debris  pieces 
obtained  in  simulation  are  all  about  0.065  kg.  Moreover,  the  counting  of  debris  in  the 
simulation  is  exhaustive.  Against  this  background,  it  will  not  be  surprised  to  see  the 
non-matching  comparison  between  the  simulation  and  field  test  data.  Nevertheless,  the 
same  trend  of  having  high  density  in  Sectors  A,  J  and  S  (in  the  direction  normal  to  the 
walls)  are  observed.  Also,  the  maximum  debris  throw  distance  agrees  well. 


(a)  Diagram  of  debris  collection  sectors 


(b)  Map  from  Kasun  data[2]  (c)  Map  obtained  from  simulation  including  twice  bounce 

Fig.  9  Comparison  of  debris  dispersion  map 

From  the  dispersion  map,  the  hazard  zone  can  be  deduced.  The  map  of  hazard  zone  is 
shown  in  Fig.  10.  The  hazard  zone  with  high  hit  probability  is  normal  to  the  walls,  and 
this  is  consistent  with  the  dispersion  map.  Table  2  gives  the  comparison  of  IBD.  The 
IBD  value  along  sector  J  is  close  to  the  field  test  result,  but  not  that  along  sector  S 
where  the  IBD  value  is  less  than  that  of  test  result.  It  is  worth  noting  that  the  field  test 
results  are  derived  from  the  assumption  which  regards  all  collected  debris  as  hazardous 
and  lethal  due  to  lacking  of  velocity  data.  This  would  lead  to  higher  estimation  of  IBD. 
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Hazardous  Zone  (Box  Shape) 
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Range (m) 


Fig.  10  Map  of  hazard  zone 


Table  2  Comparison  of  IBD  (m) 


Sector  A  (Front  wall) 

Sector  J  (Side  wall) 

Sector  S  (Back  wall) 

Global  (Max) 

Simulation 

140 

170 

160 

170 

Kasun  test[2] 

- 

180 

220 

220 

5.  Conclusion 

In  this  paper,  study  of  debris  throw  and  dispersion  after  break-up  of  reinforced  concrete 
structure  under  internal  explosion  is  presented.  The  debris  formation  algorithm,  debris 
motion  trajectory  including  the  effect  of  air  drag,  and  debris  bounce  are  investigated. 
The  debris  launch  velocity  and  launch  angle  distribution,  debris  dispersion  map  and 
IBD  determination  are  analyzed  and  compared  with  Kasun  field  test  results.  In  general, 
the  simulation  results  are  in  good  accord  with  the  field  test  results. 
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Appendix  A:  Debris  formation  algorithm 

A.l  Element  displacement  and  velocity  after  break-up 


The  information  on  displacement  and  velocity  of  all  nodes  during  break-up  process  is 
stored  in  the  output  file  named  NODOUT  in  sequence  of  intervals  Ats.  Thus  for  a  given 
time  t,  the  displacement  uECl  and  velocity  vECt  of  an  element  center  is  the  average  value  of 
its  8  corner  nodes,  which  are  represented  as 


u 


EC  t 


(A.1) 


EC/ 


1  8 


str 


(A.2) 


where  ul  and  Vj  are  the  displacement  and  velocity  of  the  ith  node  in  this  element  at  time  t. 


A.2  Debris  definition  from  the  recovered  eroded-  element  pool 

From  the  information  recorded  in  the  file  named  MSSG  for  the  eroded  elements,  the 
element  number  and  its  erosion  time  can  be  retrieved.  The  index  set  of  erosion  elements 
are  represented  by  @  .  All  the  eroded  elements  in  the  record  will  be  grouped  according 
to  the  sequence  of  time  intervals.  GEE  and  tfE  are,  respectively,  the  jth  eroded  element 
group  and  its  erosion  time. 


(j)EE  =  j/|Ele  i  has  eroded,  i=l,  2,  •••,  /VTEJ  (A.3) 

GEE  =\i\tEE -tf  <A tGE,iefE  and i  =  1,2,  ••• ,  /VTE }  (A.4) 

t?E  =  min  [tfE )  I  i  e  $>EE  and  i  =  1,  2,  •  •  • ,  NTE 

V  (A.5) 

tJE  =  min (tEE ) |zfE  > tJE  +  A tGE , i e  (j)EE  and  i  =  1,2,  •  •  • ,  NTE 

where  tEE  is  the  erosion  time  of  element  i,  AtGE  is  a  specified  interval. 

For  each  eroded  element  group,  we  can  obtain  the  number  and  size  of  debris  through 
checking  initial  inter-element  connection.  Consequently,  debris  can  be  defined  from  the 
eroded-element  pool.  It  should  be  noted  that  the  recorded  erosion  time  of  one  element 
actually  may  not  be  the  same  as  its  physical  break-up  time. 

A.3  Debris  definition  from  the  un-eroded-element  pool 

For  un-eroded  elements,  the  current  distance  of  any  nodal  pair,  which  is  at  the  same 
location  but  tied  over  the  interface  between  two  neighboring  elements  in  the  initial  state, 
can  be  obtained  from  the  results  of  last  time  step.  If  the  distance  is  greater  than  a 
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specified  maximum  distanced^  ,  the  pair  of  nodes  splits  (i.e.,  the  constrained  nodes  tie 

between  the  two  nodes  is  released).  Note  that  two  neighboring  elements  may  still  tie 
together  if  only  one  of  their  nodal  pair  splits  (or  separates)  since  there  are  4  nodal  pairs 
between  them  in  3D  (i.e.,  the  whole  interface  between  the  two  elements  does  not 
fracture).  If  fracture  does  occur,  the  two  neighboring  elements  separate.  Thereby  for  a 
given  element,  we  can  loop  through  its  surfaces  to  obtain  an  initial  fragment.  Upon 
repeatedly  looping  and  checking  the  current  disconnection  of  all  tied  nodal  pairs  along 
the  outer  surface  of  a  sub-region,  the  possible  formation  of  a  piece  of  debris  can  be 
identified,  and  eventually  a  final  form  of  fragment  (or  debris)  is  obtained  when  no  other 
elements  are  tied  to  it.  By  treating  the  un-eroded  elements  in  this  way,  we  can  define  the 
sizes  of  fragments  (or  debris),  formed  from  the  mass  of  un-eroded  elements. 

A.4  Debris  ejection  position  and  ejection  velocity 

Once  the  sizes  of  debris  are  defined,  other  physical  and  kinetic  properties  of  debris  can  be 
derived  from  the  information  of  its  component  elements.  An  individual  debris’  initial 
center  PFC( x,  y,  z),  mass  Mf,  volume  VF,  displacement  uFCt  and  velocity  can  be  defined  as 
follow. 


PK(i,y,z)  =  -f;r'£P*c(x,y,z) 
n  i= i 

nm 

i- 1 

MF=2>W 


(=1 


(A.6) 

(A.7) 

(A.8) 


in  which,  PFC(x,  y,  z),  VF,  and  P  F  ,  respectively,  represents  the  initial  center,  volume  and 
medium  density  of  the  ith  element  in  the  fragment  including  nFE  elements.  For  3D  solid 
brick  elements  used  in  this  computation,  the  initial  center  of  any  element  is  given  as 

PEC{x,y,z)  =  ^P,{x,y,z) 

8  '='  (A.9) 

where  Pfx,  y,  z)  is  the  coordinates  of  the  ith  node  in  the  element. 


Considering  a  single  debris,  it  is  reasonable  to  evaluate  its  value  at  the  centre  by  taking 
the  arithmetical  mean  of  values  from  all  its  component  elements.  The  values  of  central 
displacement  uFCt  and  velocity  vFCt  are  defined  as  follow. 


u 


FC/ 


1 


FE 


i= 1 


(A.  10) 
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V 


FCr 


1 


FE 


n 


i= 1 


(A.  11) 


By  summing  up  the  initial  position  and  the  displacements  uFCt  in  all  time  steps,  from  the 
initial  till  the  last  time  step,  one  can  get  the  final  position  of  an  individual  fragment,  and  it 
is  its  ejection  position.  Similarly,  one  can  get  the  ejection  velocity  for  an  individual 
fragment.  These  two  values  are  taken  as  the  initial  condition  in  debris  free-flying  motion. 
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Assumptions 

•  Each  debris  flies  independently 
as  a  singular  body. 

•  Ignore  spinning  motion  of  debris. 

•  Only  drag  forces  and 
gravitational  forces  acting  on  the 
debris. 


Drag  force 

F,=^CjPJV2 


Equations  of  Motion 
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V02  sin  2 60  -  V  sin  6 
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dt 
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s  =  &(F02  sin  00+V2  sin  6) 
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•  Based  on  existing  literature. 

•  Mostly  smooth  surface  with  regular 
shapes. 

•  Reynolds  Number,  range  of  interest : 
7.3  xIO3  -  7.3  xIO6 

•  Smooth  cylinder  Cd  =  0.7  (aspect  ratio 
of  2) 

•  Square  cylinder  Cd  =  2.0  (any  aspect 
ratio) 

•  Debris  with  irregular  shapes  in 
between  the  2  values. 

-^Assume  1.2<CH<2.0 
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Debris  Trajectories  by  various  Cd  (m  =  0.2  kg,  V0=  100  m/s,  Q0  =  5°) 
Final  debris  range  is  very  sensitive  to  drag  coefficient ! 
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•  Study  the  sensitivity  of 
debris  size  &  mass  on  debris 
range. 

•  Initial  launch  conditions  :  v0  = 
60m/s,  0O  =  20° ,  Cd  =  1 .2 


•The  debris  range  is  more 
sensitive  to  the  smaller  debris 
mass  (and  size)  compare  to 
the  larger  debris  mass. 

Debris  range  vs.  debris’  mass 
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Debris  bounce 
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Debris  bounce  criterion  (based  on  Knock  et  al): 

Vf  =  0;  df  =0  for  >  13°  and  l/,  >  18  m/s 

Vf=  (1.0-0.00476#,-)  Vj  and  df  =  0.484#,-  otherwise 


Incoming  debris 
Velocity  =  V; 


Case  Study 
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Based  on  one  of  the  Kasun  II  test  No  2. 

■  Comparison  of  modeling  results  with  test 

❖  Comparison  of  velocity  distribution 

❖  Comparison  of  vertical  launch  angle  distribution 

❖  Comparison  of  debris  dispersion  map 

❖  Comparison  of  IBD  (m) 
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Field  test  object  and  calculation  model 
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•  Test  conducted  jointly  by  FOI 
and  NDEAin  2006. 

•  2m  x  2m  x  2m  RC 

•  Wall  &  roof  thickness  =  0.15m 

•  (|)1 2  rebar  c/c  both  ways 

•  Bare  charge  =  20kg 

•  Loading  density  =  2.5kg/m3 


AGM-S-2.5  (KASUN-II  TEST  2) 


Air 


Numerical  model 
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Debris  percentage  (%) 


Comparison  of  velocity  distribution 
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Modeling 


Test  2 
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Launching  velocity  (m/s) 


Velocity 


Launching  velocity  by  modeling 


Launching  velocity 
by  analysis  of  video  record 
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Debris  percentage  (%) 
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Comparison  of  vertical  launch  angle  distribution 
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Modeling 


Launching  angle  by  modeling 


Test  2 
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Launch  angle 

Launching  angle 
by  analysis  of  video  record 
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Diagram  of  debris  collecting  sectors 
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Comparison  of  debris  dispersion  map 
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Test  2 
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Kasun  Test 


Modeling  with  bounce  twice 


*  Collected  debris  from  mass  bin  1 
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Debris  IBD 
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IBD:  Radius  of  the  smallesi 
circle  that  encompasses  the 
hazard  zone  totally 
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Future  works 
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•  More  test  could  be  simulated  especially  of  higher  loading  density. 

•  Investigation  works  the  debris  trajectory  and  drag  coefficient. 

Debris  cluster  effect 
Debris  spinning  effect 

•  Study  of  post  impact  problem  :  bounce,  roll  and  further  breakup  upon  impact. 

•  Effects  of  soil  cover  on  internal  loading  and  debris  throw  distance. 
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